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lation patterns and experimentally feasible BEC systems in which these effects might be detected 
in the near future. 
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I. INTRODUCTION 

Within the framework of quantum field theory in 
curved spacetime where fields propagate on a fixed back- 
ground, some interesting - and sometimes surprising - 
quantum instabilities are predicted. Among them a cru- 
cial role has been played by the so called cosmological 
particle creation: a quantum emission of particles as a 
consequence of an expansion of spacetime [TJ |5J. Be- 
ing associated with wave propagation over a dynamical 
space and not (when backreaction is neglected) to a prop- 
erty of the Einstein equations, it is not peculiar to the 
gravitational framework but appears as a quantum field 
non-equilibrium process in a diverse range of physical 
systems. This is a feature shared with other semiclas- 
sical emissions, for example Hawking radiation [3J. For 
this reason, recently there has been a growing interest 
in studying such non-trivial quantum effects in a labo- 
ratory setting with condensed matter systems, a field of 
research that goes under the general name of "analogue 
models of gravity" (see e.g. [1] for an extensive review). 

In particular, the propagation of linearized acoustic 
perturbations on an inviscid and irrotational fluid can be 
shown to be described by the same equation of motion 
as that which describes the propagation of a scalar field 
on a curved spacetime. Of course such ideal systems 
are not readily available in nature but it is possible to 
produce systems which, to a very good degree of approx- 
imation, do actually behave as perfect fluids. In recent 
years a special role in this sense has been played by Bosc- 
Einstein condensates (BEC), the basic observation being 
in this case that phase perturbations of the wave func- 
tion describing a weakly interacting BEC are, in a certain 
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limit, described by a massless quantum scalar field prop- 
agating on a non-trivially curved background Lorenzian 
space-time [JHS]- 

A BEC can therefore be used to simulate special grav- 
itational backgrounds and the associated particle pro- 
duction such as Hawking radiation of black holes [TUTU] , 
or the cosmological emission mentioned above that we 
are concerned with in this paper. For example, releasing 
the trapping potential of a BEC appropriately [TTrfTT] 
or modulating external parameters of the experimental 
setting [I8rl22j . it is possible to induce the necessary dy- 
namics which give rise to interesting background geome- 
tries. Due to the non-trivial time evolution, the system 
undergoes a non-equilibrium phase to which are associ- 
ated quantum emissions whose spectrum strictly depends 
on the parameters modulating the evolution. This opens 
up the possibility of experimentally accessing new kinds 
of phenomena directly involving the quantum nature of 
cold atom systems as well as providing additional insights 
on the general mechanisms of such emissions in the grav- 
itational case. 

Cast in a different light, this possibility of describ- 
ing terrestrial experiments using a curved spacetime ana- 
logue provides one with a new set of tools to tackle dif- 
ficult problems in condensed matter physics and is able 
to suggest new and interesting non-equilibrium quantum 
effects. Moreover, a successful application of these al- 
ternative techniques and consequent observation of the 
effects predicted would serve as further theoretical sup- 
port for the validity of the semiclassical theory. 

The issue of how to practically observe these quantum 
emissions and verify their origin is a crucial point in the 
condensed matter literature. In particular there has been 
recently a large interest in the calculation and measure- 
ment of correlation patterns in expanding condensates to 
access these quantum effects (for example see [23H25]). In 
[5] it has been proposed that evidence for analogue Hawk- 
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ing emission in a sonic black hole manifested in a BEC 
might be solicited by measuring the two point correlation 
structure associated with the quantized density or phase 
perturbation fields. The idea is that the correlation pat- 
tern should display a peculiar signature of the type of 
emission expected from black holes. Following the same 
line of reasoning, we apply this idea to analogue cosmo- 
logical emission and attempt to isolate a cosmological 
signal in the correlation pattern in an expanding BEC. 

In particular this article is concerned with the exper- 
imentally accessible correlations between the phase and 
density perturbations in a BEC which are born during 
the expansion of the underlying gas both in the widely 
studied effectively 1+1 dimensional "cigar" shaped geom- 
etry and in a homogeneous expanding BEC. Such BEC 
systems are of theoretical interest in their own right and 
correlation measurements are experimentally accessible 
with the current technology [2"3Tl25| . Therefore it is pos- 
sible that these systems will provide in the near future an 
arena within which to experimentally address interesting 
and new non-equilibrium quantum effects while at the 
same time shed light on the physics of quantum fields in 
expanding spacetime without using telescopes. 

In the first section we collect and contextualize some 
results relevant for the study of correlations in an ex- 
panding BEC, reviewing how the analogy is applied to 
BEC and the dynamics of an expanding BEC. In the 
second section we compute the correlation structure for 
a massless scalar field in a particular expanding back- 
ground spacetime in 3+1D and show how to interpret 
such a result in terms of both an isotropically expanding 
3+ ID condensate and a BEC with a time dependent scat- 
tering length. We shall develop a general formalism for 
probing the effect of particle production (Bogoliubov co- 
efficients) on correlations as embodied in the equal-time 
Wightman functions. In this context we shall derive re- 
sults which are completely generic to particle production 
in any spatially translation invariant system, which can 
be applied to standard cosmologies as much as to "ana- 
logue spacetimes" . In the third and final section which 
represents the main result of this article using both ana- 
lytic and numerical methods we discuss the structure of 
correlations in an anisotropically expanding, effectively 
1 + 1-dimensional, BEC. In particular we shall derive the 
general, finite size, 1 + 1-dimcnsional, equal time Wight- 
man function relevant for the "cigar" shaped expand- 
ing BEC presently realized in experiments. We conclude 
with some proposals for improving the observability in 
practice of the correlation signal. 



II. BACKGROUND 

A. A BEC as an analogue spacetime 

In this article we consider a weakly interacting BEC 
trapped in an attractive harmonic potential. Given the 
by now standard usage of this system as an analogue 



model of gravity we limit ourselves to a summarry of the 
salient equations which we shall need in what follows, 
directing the reader to the extant literature for further 
reading [3HS]- 

In the dilute gas approximation, one can describe a 
Bose gas through a quantum field \& satisfying 



V 2 + K x t(x)+ 5 (a)§t§ §. (n.l) 



2 m 



Here m is the mass of the constituents, g parametrizes 
the strength of the interactions between pairs of bosons 
in the gas and a is the scattering length of the atoms in an 
s-wave approximation. The former can be re-expressed 
in terms of the scattering length as 



1(a) 



4-wah 2 



m 



(11.2) 



Adopting a mean field approximation, the quantum 
field can be separated into a macroscopic (classical) con- 
densate (called the wave function of the condensate) 
and a fluctuation: & = ip + <p, with (vp) = ip. Then, 
by adopting the self-consistent mean field approxima- 
tion and neglecting backreaction effects of the quantum 
fluctuations on the condensate, one obtains the Gross- 
Pitaevskii (GP) equation of motion for the wave function 
of the condensate 

Now, adopting the Madelung representation for the wave 
function of the condensate [2S1 



V>(i, x) = y/p(t,x) exp[-i0(t, x)/ft], (II.4) 

where p = ip*tp is the condensate density, and defining 
an irrotational "velocity field" by v := V0/m, the Gross- 
Pitaevskii equation can be rewritten as a continuity equa- 
tion plus an Euler equation: 



d_ 
dt 



TO aT + v v 2 



p + V • (pv) = 0, 
+ y cxt (i,x) + .9p- 



h 2 V^y/P 

2m J~p 



(II.5) 



(II.6) 



These equations are completely equivalent to those of an 
irrotational and inviscid fluid apart from the presence of 
the so-called quantum potential 



V, 



quantum 



= -ri 2 v 2 vV(27Vp), 



(II.7) 



which has the dimensions of an energy. It is then clear 
that when small gradients of the background density are 
involved, this term can be safely neglected and a hy- 
drodynamical approximation (also known as the Thomas 
Fermi (TF) approximation) holds. 
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We are concerned with the quantum excitations on 
top of the BEC background solution of the GP equa- 
tion (or hydrodynamic equivalents). Inserting the mean 
field ansatz in Eq. ( II. 1 ) and again neglecting backreac- 



tion, one obtains the equation for the operator valued 
perturbation 

d ( ft 2 \ 

ift—(p(t,x) = ( -— V 2 + y cxt (x)+2.g ( oJ £(t,x) 



gpcpt(t,x). 



(II.; 



The elementary excitations of this dynamical system are 
described by variables obtained by a Bogoliubov trans- 
formation from the fundamental field. This can be also 
shortcut by adopting the so called quantum acoustic rep- 
resentation [3] 



<p(t,x) =e 



1 



Pi - 1 



IP 



(II.9) 



where p\ , 9\ are real quantum fields which describe the 
excitations in density a nd ph ase respectively. Inserting 
the above ansatz in Eq. (II. 8 1 one again obtains a pair of 



hydrodynamic-like equations for the quantum excitations 
d t px + ~V ■ ( p x V6» + p Wi) = 0, (11.10) 

TO V / 

v$ ft 2 

dA + : + g(a) pi- — D 2 p x = 0, (11.11) 

to 2m 

where D 2 represents a second-order differential operator 
derived from the linearization of the quantum potential 
[3] . In a long wavelength approximation this can be safely 
neglected and one can then show (in complete analogy 
with perfect fluid systems) that the perturbations 0\ do 
satisfy a second order equation of a Klein Gordon type 
for a massless scalar field propagating on a curved back- 
ground geometry described by the metric 



ds 2 = 



v 2 )dt 2 - 2v i dx i dt + dx 2 ] . (11.12) 



Here c 2 (t,x) = gp(t,x.)/m is the squared local sound 
speed of the perturbations and v(£,x) is the background 
flow velocity of the BEC introduced above [27] . It is 
then clear that in this long wavelength regime the quasi- 
particles satisfy a relativistic, phononic dispersion rela- 
tion oj 2 — c 2 k 2 (with lj and k equal respectively to the 
energy and wavenumber of the quasi-particle) . Alterna- 
tively one can adopt an eikonal (short wavelength) ap- 
proximation and show that the complete dispersion re- 
lation is indeed the well known Bogoliubov dispersion 
given by u 2 Bog = c 2 k 2 + hk 4 /(2m). It is evident at this 
point that in this context the gravitational analogy is 
valid only in the long wavelength approximation. More 
specifically, for wavelengths of the order of the healing 
length £ = ft/^rngp, when the quartic term becomes 
comparable to the quadratic one in w| 09 (fc), a geometri- 
cal description is no longer available. 



B. Scaling solution 

Condensates undergoing expansion or contraction as 
a result of a time varying confining harmonic potential 
admit a dynamic generalization to the static TF approxi- 
mation for the GP equation known as the scaling solution 
[251 129] . In this case the external potential is of harmonic 
type with time varying trapping frequencies 



1 



i;. x( = - ^muj 2 {t)x 2 



(11.13) 



Here the density undergoes a simple scaling transforma- 
tion in each spatial dimension Xi parametrized by the 
so-called scaling functions Xi(t) as 



p{t, x) = 



p s (x l /\ i {t)) 



where p s is the static TF ground state solution^] 



(11.14) 



Po 1 



if ^CoiX 2 <\, (11.15) 



p s = otherwise, 



(11.16) 



where LUi — muj 2 /2p (/i being the chemical potential) and 
Po is the maximum density at the centre of the BEC. 

In the static TF approximation the density and chem- 
ical potential can be expressed algebraically in terms of 
the physical parameters as 



Po 



p _ 1 fl5m 3 NU t uj^ 2/5 



g(d) 8tt V a 3/2 ^ 3 



and 



1. / _ AT _ mw 
p 2 V V ft 



2/5 



(11.17) 



(11.18) 



respectively, where N is the number of condensed bosons. 
Here we have introduced the notation w — (Yli^i) 1 ^ i n 
d spatial dimensions, the geometric mean of the trapping 
frequencies. For the phase, the scaling transformation 
takes the form 



9(t,x) — -to x 2 \, \j 



(11.19) 



where ' represents the lab time derivative dt ■ 

This scaling induces a background flow velocity profile 
in the condensate given by 



v(i,x) 



AM 
X(t) 



(11.20) 



1 Note that the TF approximation is not a good approximation 
near the boundary of the condensate where the density drops 
to zero. The condensate boundary, here given by the parabola 
sharply dropping to zero, in the full solution will be smooth. 
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With this velocity and the evolution of the density given 
in (11.14) one can start studying the analogue geometry 



associated with s uch a scaling background solution using 
the metric ( |II.12| ). 

Of particular experimental interest is the "cigar" 
shaped condensate possessing one elongated dimension 
("||") and two tightly trapped orthogonal dimensions 
("_L"). In what follows, unless otherwise stated, we work 
with the cigar geometry. Since our main focus will be 
on such elongated condensates, we list here the following 
realistic values for the physical parameters describing the 
cigar BEC 



m - 
a 

X : 



= 10 2 Hz, w_L = 10 3 Hz 
= 1.57 x l(T 25 kg (Rubidium 87) 
42 Bohr = 2.221 x 1CT 9 m 
= 10 5 . 



(11.21) 
(11.22) 
(11.23) 
(11.24) 



These are respectively the trapping frequencies, atomic 
mass, scattering length and number of condensed atoms. 
Assuming these numbers, the following derived quantities 
are obtained 



/iR!5.8x 10 
po ps 2.9 x 10 



-31 



20 m -3 



5.45 x 10~ 5 m (11.25) 
kgmV 2 fa 5500/i s" 1 (11.26) 
(11.27) 

c S:0 ~ 2 mms -1 (11.28) 

corresponding to the transverse condensate size (where, 
in the TF approximation the density drops to zero), 
the (static) chemical potential, the maximum conden- 
sate density and the maximum initial sound speed re- 
spectively. 

For time dependent harmonic traps the scaling func- 
tions Xi(t) are not free parameters but are constrained to 
be solutions to the auxiliary equations of motion 



Ai(t)+w 2 (t)Ai(i) 



^ 2 (0) 



(11.29) 



with boundary conditions Ai(0) = 1 and Ai(0) = 0, de- 
termined by a self consistency condition on the dynamics 
of the expanding condensate [25] ■ 

Starting from such a cigar shaped initial geometry, let 
us consider the case when the "||" potential is suddenly 
switched off while keeping the "_L" potentials constant. 
A similar problem has been studied in [53] and observed 
in [23] where the authors instead consider a BEC freely 
expanding in all the three spatial dimensions after com- 
pletely releasing the trapping potential of an initially 
cigar shaped condensate. In these analyses, due to the 
expansion in all three dimensions, the density drops very 
quickly and the gas becomes collisionless signalling the 
failure of the hydrodynamical description shortly after 
the release of the trap. In the case we consider here in- 
stead, the perpendicular trapping is kept tight and the 
density drops more slowly and the condensate remains 
longer in a collisional regime where the spacetime anal- 
ogy pertains. 
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FIG. 1. The solution of the scaling problem for an anisotropic 
expansion of an initially cigar shaped configuration for both 
the radial/orthogonal (blue-dashed line) and axial/parallel 
(red-solid line) directions to the axis of release. Note that 
the scaling function of the direction in which the trapping 
frequency remains fixed (blue, dashed) in fact decreases dur- 
ing the expansion. The parameters used in this solution are 
u\\ — 100Hz, uj± = 10 3 Hz and a sudden release of the || di- 
rection only. 



In Fig. [T] we show the numerical solution to (11.29) for 
cigar- like initial geometry after the potential in the "||" 
direction is suddenly switched off while keeping the "_L" 
potentials constant. Notice that the tightly constrained 
"_L" dimension in fact pulls even tighter during the ex- 
pansion in the "||" dimension as indicated by the decay 
of the scaling functions Aj_. In this way one can under- 
stand how the volume of the cigar geometry, scaling as 



(A||Ai 



decays much more slowly than in an isotropi- 



cally expanding BEC which scales as An 3 . 



C. Validity of the hydrodynamic approximation 

As already mentioned above, it is also important 
to keep track of the hydrodynamical approximation 
throughout the expansion. Since the quantum pressure 
scales differently with time than the density, it is of inter- 
est to understand whether the hydrodynamic approxima- 
tion improves or degrades as the expansion proceeds and, 
if it loses accuracy, when the approximation breaks down. 
In Fig. [2] we display the quantum pressure V^uantum and 
interaction VS n t er action (= 9p) contributions to CP before 
the release of the axial trapping potential. In Fig. [3j we 
display the ratio V'qu^tum/lmteraction as a function of co- 
moving axial position z = zAm for four different times 
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Quantum pressure vs interaction term at t=0 
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FIG. 2. The quantum pressure and interaction energy ex- 
pressed in Joules as a function of axial distance in meters. 
Note that the quantum pressure is certainly non-negligible 
near the boundary of the TF approximate ground state den- 
sity wavefunction (at z ~ 2.7 x 10 _5 m). 
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FIG. 3. The ratio R = Vquantum/Vinteraction as a function 
of co-moving axial position 5 = zAn at various times t after 
release of the trapping potential. Note that the axial half- 
length of the condensate at time t = (which of course is the 
co-moving length for all times) is L z ~ 2.7 x 10~ 5 m. 



after the expansion. 

Note that although there is always a region near the 
edge of the condensate where the ratio is greater then 
1, as time progresses the proportion of the condensate 
which remains within the hydrodynamical regime gets 
smaller until approximately t — 0.25s when the quan- 
tum pressure becomes important globally in the conden- 
sate. The failure of the hydrodynamic approximation 



near the boundary of the condensate is expected in gen- 
eral since the TF approximation is known to break down 
on scales on the order of a healing length from the bound- 
ary where the density becomes exponentially small. In- 
deed the true exact solution, being a smooth exponential 
decay to zero, differs from the TF approximate solution, 
given in (11.16 1, near the edge which sharply goes to zero 
there. Our analysis in the following will be concentrated 
on the more central part of the condensate where the TF 
approximation holds. 



III. QUANTUM EXCITATIONS 

We now intend to consider particle creation associated 
with the expansion of the condensate and its signature 
in the correlation structure. Before doing so however, 
one might wonder if such expansion is sufficiently rapid 
to lead to any relevant particle production at all. Also, 
given that we shall work in the analogue gravity frame- 
work, one would also like to be sure that the ostensible 
excitations can be meaningfully described as phonons 
rather than having to deal with the complicated issues 
associated with non-linear dispersion. We shall discuss 
these issues first and then introduce the correlation func- 
tion that will be the subject of our studies. 



A. Timescales 

Let us now analyze and compare the relevant inverse 
time scales involved in the problem of an expanding 
condensate to make an order of magnitude estimate of 
whether excitation effects can be expected to be observ- 
able. 

Firstly there is an intrinsic infrared cutoff given by the 
size of the condensate. This length scale translates into a 
frequency u s using the (central maximum) sound speed 
u s = c s /(L\(t)). 

Secondly there exists an intrinsic ultraviolet cutoff 
associated with the healing length, the healing fre- 
quency Wh = 27tc s /£ where £ = h/(mc s ) is the healing 
length. The healing length, as already mentioned, is the 
length scale at which the quantum potential, neglected 
in the hydrodynamical approximation, becomes intrinsi- 
cally comparable with the interaction energy (which is 
proportional to the density of the condensate and hence 
becomes less relevant at lower densities). This is also 
the scale at which the Bogoliubov dispersion relation en- 
ters in the quadratic regime and hence the description in 
terms of non-interacting phonons breaks down. 

Thirdly there is the characteristic frequency associated 
with the expansion itself, w e ~ A/ A. which fixes the fre- 
quency of the most abundantly created excitations in an 
expanding condensate. This is generic to any particle 
creation effect by a time-varying external field (a "dy- 
namical Casimir effect") see e.g. |18) . 
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It is of interest to compare these three frequencies in 
order to test the validity of our approximations and to 
estimate if the expected excitations are created within 
the natural ultraviolet and infrared cutoffs. In Fig. [4] we 
compare these three inverse time scales using the accu- 
rate parameters listed above. 



Expansion and size frequencies 
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FIG. 4. A comparison of the relative magnitudes of uj e and lj 3 
(red-solid and blue-dashed respectively) and with the healing 
frequency (inset plot) in s" 1 as a functions of time (in s). Note 
that, for the majority of the duration of the expansion, ui s is 
much smaller than the expansion frequency cj e which is much 
smaller than the healing frequency. Hence the bulk of the 
particle production spectrum lies within the low frequency 
cutoff provided by the total size of the condensate and the 
high frequency cutoff provided by the healing length. 

We see that indeed there exists a regime in which the 
typical excitation frequency lies entirely within the natu- 
ral ultraviolet and infrared cutoff frequencies. Quantita- 
tively this regime begins after about 0.005s. Let us also 
note the vastly higher healing frequency representing the 
UV cutoff for a phononic treatment. This is shown in the 
smaller subhgure in Fig. [4] for comparison. 

Therefore we conclude that there is a viable regime 
within which the excitations have a phononic (linear) dis- 
persion relation and wavelengths much shorter than the 
size of the BEC. 



B. Wightman functions 

The Wightman function G(x, y) := ((j)(x)<f)(y)) for a 
quantum field (/> is a measure of the correlation between 
two separate points x and y in spacetime. If the field 
<j) were a truly random variable the correlation would 
be zero everywhere. In general the equal-time Wight- 
man function for a massless scalar field on a flat 3+1 
dimensional FRW spacetime with scale factor a (which 



includes Minkowski space as the a = 1 special case^ 2 ] 
only depends on the magnitude of the co-moving spatial 
difference x — ||x — y|| and (possibly) on time. This ho- 
mogeneous function thus is really only a function of these 
two variables 

f°° dk k 
= ^ - sin (fcr) u k (v)u* k (v) 
Jo zi7r x 

= G(r,,x) (III.l) 

where k = Vk • k, it& are the mode-function solu- 
tions to u'l + (fc 2 — a"/a)uk = and = (frfeUfc + 
$_ k u\) I \\[2a(r\)) , prime denoting the conformal time 
derivative d v . In Minkowski spacetime the integral is 
computed exactly as G{q,x) = H/{^-k 2 x 2 ) for non-zero x 
while on the light-cone the Wightman function also pos- 
sesses a singular imaginary part [30] . In 1+1 dimensions 
the equal time Wightman function becomes logarithmic 
and time independent G(x,y) — —h/4ir In (x — y) 2 and, 
owing to the conformal invariance in 1+1 dimensions for 
the massless scalar field and the conformal flatness of 
FRW spacetimes, also takes this form in an arbitrary 1+1 
FRW background in terms of co-moving coordinates and 
the conformal time variable [2]. 



IV. CORRELATIONS IN A ISOTROPICALLY 
EXPANDING BEC 

The main purpose of this section is to gain an intuitive 
understanding of correlations due to particle production 
and provide a useful toolbox for the more challenging in- 
vestigation of the anisotropic scaling solution. In order 
to do so we shall study first the ideal case of a homoge- 
neous 3+1D isotropically expanding BEC. In this case, 
the analogue metric is that of a simple FRW spacetime 
geometry. We shall choose a scale factor that increases 
from an initial constant value to a final constant value 
a f smoothly as in this case the structure of correlations 
can be solved for analytically. Our interest in these ex- 
actly soluble models derives from the fact that in these 
systems one can see the appearance in the Wightman 
function of the characteristic features expected from the 
creation of particles in a transparent way. Finally, we 
shall also discuss how to simply interpret such models in 
terms of a BEC with a varying scattering length. 



2 Recalf that a FRW spacetime is described by the metric ds 2 = 
a 2 {rj) {—dry 2 + dx 2 ) here written in co-moving coordinates x and 
conformal time rj where a(rf) is a scalar function of time. 
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A. Exactly Soluble Expanding Spacetime 



Written in terms of co-moving coordinates 



Xi/Xi(t) the metric (11.12) for the scaling solution BEC 
reads 



ds 2 



'p(*,x) 
gm 



(-c 2 s (t,x)dt 2 + \}{t)dZ l2 ) (IV. 1) 



from which a number of standard representations for a 
FRW background are easily obtained. The cosmologi- 
cal assumption of homogeneity is satisfied for a BEC in 
the central region where the density is almost constant. 
Working in such a central region, p and c s become func- 
tions only of t and the metric can be written in standard 
FRW form as 



ds 2 = -dr 2 + a 2 (r)dx 2 , r(t) = J dt 
with the scale factor defined by 



m 3/4 

(IV.2) 



a 2 {T) = M±\ 2 {T)^^W) 
V 9 m 



(IV.3) 



where we have used the scaling solution p(rj) = po/X 3 (rj) 
for the density. 

In what follows however, we will work with the form 
for the metric given by 



ds 2 



-a 6 (r]) drj 2 + a 2 {r ] )dx 2 



which can be obtained form the form (IV.2) by the coor- 
dinate transformation 



/ m \ 3 / 2 i 



(IV.5) 



The functional relationship between the scale factor a the 
scaling function A is unchanged in this representation as 
a 2 (77) oc y^Xjq). 

Furthermore we will consider the specific form for the 
scale factor 



4 1 4 4 4 
a i + a f , a f ~ a i 



tanh I — 
2 Wo 



(IV.6) 



where 770, cij and a/ are constant numerical parameters. 
Such a metric describes an isotropic homogeneous space- 
time with flat spatial sections which expands by a finite 
amount over a finite time period. This spacetime has 
been studied previously in the literature for example in 
[18) where it was shown, to arise as the analogue geom- 
etry experienced by phase perturbations in a BEC with 
time varying interaction strength g. In that article the 
particle content of the initial Minkowski vacuum state 
after expansion is studied in detail and the observabil- 
ity of the created particles discussed. This spacetime is 



also studied in [16] where it is shown to arise as the ana- 
logue geometry associated with a ring of trapped ions. 
A very similar geometry is studied in [52] where instead 
of tanh, the authors employ an error function^] temporal 
profile for the variable scattering length and the correla- 
tion structure is studied semi-analytically in the so-called 
"sudden limit" (which we will discuss later) , and numer- 
ically in the general case, with the full Bogoliubov spec- 
trum. In [31] the author studies the controlled release 
of a condensate such that the trapping frequency varies 
exponentially in time and the associated particle produc- 
tion in terms of the number of particles produced. The 
separate case of an exponential variation of a time de- 
pendent coupling constant g(t) is also studied as well as 
the mixed analysis where both the trapping frequency 
and coupling constant are varied again in terms of the 
number of particles produced. 

The field equation for a massless minimally coupled 
scalar field <f> propagating on the background geometry 



(IV.4) decouples into the independent mode equations 



a 4 k 2 4>k 







(IV.7) 



where 4>k — (2tt) 3 ^ 2 J d 3 ke lk ' x is the kth Fourier com- 
ponent of the field <p. For a{rf) as (IV.6) the field equa- 



tion (IV.7) is solved by a rather complicated product of 
elementary and hypergeometric functions [18j . Conve- 
niently for us, it is not necessary to work explicitly with 
these cumbersome functions: For each co-moving mo- 
mentum k the physics of cosmological particle produc- 
tion is contained only in inner products (which we label 
(IV.4) as ak and (3k) between two particular solutions <p k n and 



k to (IV.7), solutions which converge respectively in 



the past (77 —J- —00) and future (77 — > 00) to plane waves^] 
(see for example the standard monograph [2]). 

Physically, these two particular solutions represent a 
description of two different vacuum states, |in) and |out), 
in the two asymptotic regimes where the function 0(77) 
becomes constant and the definition of particle state is 
unambiguous. For example, the quantum state |in) de- 
scribed by the choice <f> k n (the vacuum state for 7/ — » —00) 
is no longer a vacuum (zero particle) state as rj — ¥ 00 since 
that state, |out), is described by the mode functions <j)° k ut 
and in general one has 



(IV. 



with /3fc 7^ 0. The spectrum of particle content of the 
state I in) at late times in terms of the late time parti- 
cle states is given by the modulus \fik\ 2 which has the 
characteristic bell shape shown in Fig. [5] here shown as 



3 Recall that the error function is defined as Erf(t) oc 
/' exp(— x 2 ) dx and smoothly interpolates between constant val- 
ues. 

4 Specifically afc is the innerproduct between positive frequency 
plane waves whereas 0^ is the inner product between one positive 
and one negative frequency plane wave 
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FIG. 5. The characteristic shape of the spectrum of produced 



particles in an expanding spacetime described by (IV.4I, as 
a function of the modulus k — vk ■ k. Recall that the 
spectrum is proportional to the volume of space V/2n as 
Nk = Vk' 2 pf./2-K so that what we plot here is in fact par- 
ticle number density. The peak wave number is fixed by the 
typical timescale of the expansion k pea k oc l/ryo- 



a function of the isotropic modulus k = Vk • k appro- 
priate for the isotropic mode equations (IV.7) and hence 
isotropic coefficients a k and fi k . 

The inner products a k and /3 k (known as Bogoliubov 
coefficients) are given exactly in this model by the simple 
expressions 



Oik 



2-/AB Y(—iA)Y(—iB) 
A + BT 2 (-i(A + B)/2) 

-2y/~AB T{-iA)T(iB) 
B-A T 2 {i{B-A)/2) 



(IV.9) 



where T is the Euler function and the dimensionless vari- 
ables A and B are 



A = krj a 2 , B — kr] a 2 . 



(IV.10) 



They are constrained to satisfy |afe| 2 — |/3fc| 2 = 1 by the 
constancy of the Klein Gordon norm of the (unique) so- 
lution 4> written in the two different mode function bases. 

Since the mode-functions are simple linear combina- 
tions of plane waves long after the expansion has taken 
place, the equal time Wightman function there can be 
written down in a relatively simple way in terms of a k 
and 0k alone as 



47T 

ztt / JO 



where 



dk 



sinfca; 



l + 2\/3 k \ 



sin kx 



i + r 

(IV.11) 
(IV.12) 
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FIG. 6. The exact numerical integration of (minus) the entan- 
glement contribution to G as a function of point separation 
x for four successive times rj = (red-solid), rj = 1 (green- 
dashed), 77 = 2 (blue-dotted) and r\ = 4 (black-dot-dashed) 



for the choice of parameters a i 



4 



-K and 770 = 1. 



Note again that, given spatial homogeneity and isotropy 
of the system, here G is a function only of the magnitude 
of the spatial separation. Of course, such an expression 
cannot be considered exact for a realistic finite-sized BEC 
but will break down nearby the boundaries due to finite 
size effects. Furthermore, in the case of a finite volume 
"cigar" -shaped BEC and anisotropic expansion the ex- 
pression (IV. 11 1 will have to be replaced by a more com- 
plicated formula as we shall see in sect. |V| 



The expression (IV.ll) is the sum of three physically 
distinct terms. The 1 term represents the background 
zero point fluctuations, the |/3fc| 2 term an enhancement 
of the background correlations due to particle production 
and the terms labelled / and I* represent additional non- 
trivial structure propagating on top of the (enhanced) 
Minkowski-like correlations due to the entanglement be- 
tween the created pairs of particles. The / + /* term 
can be integrated exactly as an unwieldy infinite sum of 
functions which we do not write down here. In Fig. [6] we 
plot the result of numerically integrating the contribution 
from I + r 



47T 



dk (I + I* 



(IV.13) 



in G for four successive times as a function of separation 
x = ||x — x'||. We see that going backwards in time a 
bump of correlations merges into the correlation singular- 
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itj{^] at the origin. The (coordinate) speed of the bump 
is determined by the final value of the scale factor a 2 
which in the displayed plot was taken to be n. We note 
that the bump moves with a velocity of close to 7.5 units 
per second which is approximately the speed at which 
two particles moving at speed tt units per second would 
be receding from one another. Note that the coordinate 
speed is not constrained to be the speed of sound c s , the 
coordinate speed in lab variables. 

In the next section we show that the essential structure 
of these correlations are present already in the so-called 
adiabatic approximation to the integrand I + I* . 



1. Adiabatic and sudden approximations 

The adiabatic approximation for a mode k is accu- 
rate when krjQ 3> 1. It is a good approximation for 
high wave number, or short wavelength modes. These 
are the modes which feel the expansion as a "slow" pro- 
cess. In the adiabatic limit one can show, from (IV.9) 
and pVT2|, that 



I+I* ~ - sin (kx) sin [k ($ - 2a 2 f T])] e^™^ (IV.14) 
where the phase is given by 



$ = In 



In 



V 



a 2 f r) - (IV.15) 



On the other hand, the sudden approximation is accu- 
rate for a given k whenever krjo <C 1. Hence it should be a 
good approximation to the exact result for low wavenum- 
bers. In this sudden limit one can similarly show (again 
from and ||IV.12j» that 



I + I* 



la) 



a>i sin kx 



2 2 

a 2 a 2 



cos(2fccif77) 



(IV.16) 



which represents an unbounded contribution to G if as- 
sumed accurate over the entire positive k axis (as would 
be the case in the formal limit ?/ — >• 0) . In Fig. [7] we 
compare the exact, adiabatic approximate and sudden 
approximate forms of I + I* for a representative choice 
of x and r\. We note that the adiabatic approximation 
is accurate almost uniformly over all k, even close to the 
origin, where it is superseded by the sudden approxima- 
tion while still capturing the qualitative behavior of the 
I + I* term. 



5 Since we always work under the assumption that the field the- 
ory is an effective field theory valid only above a small length 
scale, for example the healing length in the case of a BEC, the 
coincidence singularity does not exist in practice. In the BEC 
example, higher order contributions to the dispersion relation 
become non-negligible and in fact dominate at very short length 
scales resolving the singularity. 
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FIG. 7. A comparison of the exact expression (in red-solid) 
with both the adiabatic approximation (in blue-dots) and the 
sudden approximation (in green-dashes) near the turn over 
region (kr/o w 0.3) where neither approximation is very accu- 
rate. We use the representative values of t = 2 and x = 3. 
Note that the sudden approximation is most accurate for 
small k and the adiabatic for large k. Also note that the 
adiabatic approximation remains relatively accurate (at least 
in mimicking the qualitative behavior of the exact function) 
all the way down to k = where the approximation in prin- 
ciple should be poor. 



Integrating the adiabatic approximation I + I* over k 
we observe a "bump" of correlations traveling out from 
the origin and decaying as displayed in Fig. [8] This is in 
close agreement with the exact result as shown in Fig. [6] 
We note that the only feature of the exact result not 
captured by the adiabatic approximation is the negative 
"tail" to the bump as it propagates outwards from the 
origin. 

The role of the parameter 770 is visible in Fig. [9] where 
we again show the result of integrating the adiabatic ap- 
proximation^] to I + 1* for the same values a 2 — 1 amid 



= tt and times rj = 0, 2 , 4. One sees that as the rate 
of expansion increases (as 770 becomes small) the shape of 
the bump becomes sharply peaked and that the peak oc- 
curs slightly closer to the origin. This is consistent with 
our intuitive understanding of the bump as arising from 
propagating particles created during the expansion. 

For small 770 the expansion and hence the particle pro- 
duction takes place over a short time interval. Assume for 
the purposes of illustration that all the particle produc- 
tion occurs at a one point in space, say at X. (Particle 
production occurs at every point but this simplification 
serves to capture the relevant physical mechanism) Then, 



For computational reasons we restrict our attention to the adia- 
batic approximation in this section. 
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FIG. 8. The integrated I + 1* contribution to the Wightman 
function in the adiabatic approximation at three successive 
times r\ = (red-solid), rj — 2 (blue-dots) and r\ — 4 (green- 
dashes). Note that since the expression ( IV. 1 1 1 is only valid 



after the expansion has taken place (the midpoint of which 
occurs around r\ = 0) the contribution for t < are not 
shown. To make these plots and for comparison with the 
exact result again we have used the values af = 1, a* = tv, 
and r/o = 1. Also plotted (in black- dash-dot) is the standard 
flat spacetime correlations (Aw 2 x 2 )~ 1 . 
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FIG. 9. Comparison of the bump structure of correlations for 
three different values of the parameter 770 = 2, 1, 1/2 and 1/4 
(from top left to bottom right) at three different times 77 = 
(red-solid), r\ = 2 (blue-dots) and 77 = 4 (green-dashes). Also 
plotted (in black- dash-dot) is the standard flat spacetime 
correlations for comparison. 



two observers at a fixed spatial separation symmetric to 
X will simultaneously observe a burst of particles pass- 
ing by. For large 770 the particle production is spread 
out over a longer time scale since the scale factor a is 
changing significantly over a longer period and hence the 
observers will see a more diffuse burst of particles. In this 
way it is not surprising that attempting to integrate the 
sudden approximate solution over all k (that is, in the 
formal 770 — > limit) one finds a propagating singularity 
at approximately x = 2a^?y as can be inferred from the 



structure of (IV. 16). 



2. Varying scattering length interpretation 

As mentioned at the beginning of this section, the met- 
ric (IV.4) can be used to describe phase perturbations 



in a BEC with a time dependent scattering length a(t) 
(and hence coupling parameter g(t) and sound speed c s ) 
as discussed in the article [TS] . In such varying scattering 
length models one works in a homogeneous approxima- 
tion where the background condensate density p is as- 
sumed constant and the flow velocity v everywhere zero. 
The FRW scale factor (which we temporarily write as a s 
for clarity) is related to the scattering length a by 



aM = (4^) 



1/4 



(IV.17) 



Hence a decrease in the scattering length corresponds to 
an expansion of the FRW spacetime and vice versa. In 
ITSl the authors use this form for the metric and also the 



tanh scale factor profile (IV. 6 1 described above. 



We interpret our result in this context as a propagating 
bump of correlations when the scattering length is varied 
according to 



a(ry) 
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4ttH 2 



tanh 




(IV.18) 

Note that this time dependence is not determined by any 
internal dynamics of the BEC but instead is a free ex- 
perimental choice. This is in contrast with the expanding 
BEC case where the time dependent parameter, the scal- 
ing function A, is determined by the internal dynamics of 
the BEC. 



V. AXIAL CORRELATIONS IN AN 
EXPANDING ELONGATED CONDENSATE 

We shall now consider the experimentally relevant case 
of an anisotropic harmonically trapped BEC, where the 
trap is released only along the z dimension. 

We write A(t) := \ z (t) and Aj_(t) := A x (t) = 
X y (t). The simple form and cylindrical symmetry of 
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the scaling solution in this case allows us to separate and 
the variables in the background velocity as c s (t, z, r) = 
c{z,r)/y / \(t)\ 2 L (t) wher e c = c s Vl - u z z 2 - uj r r 2 and 6Q(t,z) 
u)i was defined in section [11 B| 

Let us now introduce new coordinates through the ex- 
act differential expressions 



dT 



dt 



A 3 / 2 A ± ' 



dZ 



dz 
c(z, f) 



fdr 



with f(f, z) = — J z d?c/c 2 dz' . Then the metric (IV. 1 1 is 
written as 



ds 2 = 




~ dT 2 + dZ 2 



A5 



de 2 




2/ dZdr 



(V.2) 



Note that / vanishes in the limit where we can neglect 
the radial derivatives of c. In that limit note that the 
metric takes a particularly simple diagonal form. 

Relabeling <f> :— B\ the action for phase perturbations 
written in these coordinates is 



S=- 



i J d^^g^d^O^ 
~ J dTdZ [-n(T,Z)(d T cbj 2 + H(T,Z)(d z <j>) 2 ] 

+<D(d r (f>) + 0(dgdp) (V.3) 



where 



and 



n = / drde 




H= f drd6.p^r\\ 

J \ gm 



/Af ~ 
A i 



(V.4) 



n + sn. (v.5) 



The action (V.3) has been written up to terms involving 
the radial and angular derivatives of the field <fr which 
we will discard at this stage. The justification for this 
approximation is transparent: We consider condensates 
sufficiently tightly trapped in the radial direction such 
that the infrared cutoff provided by the finite radial size 
and the ultraviolet cutoff provided by the minimal heal- 
ing length in fact preclude the appearance of any radial or 
angular modes at all. In other words, in order for a mode 
of the field <\> to satisfy the ultraviolet band-limitation it 
must necessarily not contain any non-zero radial or an- 
gular frequencies. Or even more simply, due to the small 
radial size of the condensate radial and angular modes 
necessarily would oscillate on length scales shorter than 
the healing length and are hence absent. 

Explicit formulas for the functions fl and SQ in terms 
of the lab coordinates (t, z) are obtained by integrating 
(TO) and dvT5b 



Q(t,z) 



fio Aj 



Vx 



i 



A 2 



3/2 



(V.6) 



2» A 3 / 2 
ujj_ Aj_ 



--1/2% Q Z 

"ii x~ 3 Y 2 



-3Q 



1/2 Z 

II A 3 



(V.l) where 



fir 



-in 
3ff 



(-) 



3/2 



(V.7) 



(V.8) 



Here it is understood that these functions are zero if the 
right hand sides become negative. Both functions f2 and 
SQ have support only over the extent of the condensate. 
n monotonically decreases from its maximum at z = 
whereas <5f2 vanishes at z = increases to a maxima at 



--1/2 

' J \\ 



/4 before decaying back down to zero. With 
the parameters we are using here the maxima of 5Q is 
approximately £lo/(u>iiU)±) which, due to the large trap- 
ping frequencies taj, is five orders of magnitude smaller 
than CIq/u!j_, the maximum of Q. Hence in what follows 
we set H = il, neglecting the <5f2 contribution to H. 

In this limit (and neglecting terms containing radial 
and angular derivatives, as discussed) the 1 + 1 dimen- 
sional action (V.3) becomes particularly simple 



S=- 



X - J dTdZVl(T, Z) [-(Or, 



+ {dz4>f], (V.9) 



leading to the field equation 

- d T (Qdr^) + d z {Q,d z <t>) = . (V.10) 
In what follows we shall work with this field equation. 

A. Conformal symmetry approximation 

As a zeroth order attempt at uncovering a signature of 
"cosmological" particle creation in the Wightman func- 
tion one might be tempted to discard derivatives of the 
function Vt. Under such an approximation the field equa- 
tion reduces to the flat 1+1 dimensional wave equation 
for the scaled field if) = Sl<j> whose solutions are simple 
exponentials. The Wightman function for the phase field 
<f> is then simply related to the standard Green function 
as 



(<MT,Z)0(T',Z')) = 



n(T, Z)fl(T',Z') 
h In [Ax + Ax~] 
Q(T, Z)fl(T\ Z') 



CV.H) 



where 



are the characteristic null coordinates. How- 



ever this is immediately seen to be too strong an approxi- 
mation to capture any particle production effects. In fact 
the above correlator can be easily plotted as in Fig. [TO} 
which shows no sign of the typical transient propagating 
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Field correlator during eternal expansion as a function of 
time 
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FIG. 10. The field (phase) correlations G"* as a function of 
time t for three different fixed laboratory spatial separations 
0.2L Z , 0.5L 2 /2and 0.8L x /2. 



over-correlation feature normally associated with particle 



production (see section IV A) 



This observation seems to be at odds with several re- 
sults found in some the literature concerning the corre- 
lation structure in the case of Hawking radiation. For 
example, in [5J [S] an analogous calculation is performed 
to the one here, dimensionally reducing a 3+1D BEC 
dynamical problem down to one in 1+1D in the case of 
an acoustic black hole background geometry for phase 
perturbations. The important point for the current dis- 
cussion is that the authors of [El [9] in fact do observe 
numerically a signature of created Hawking quanta in 
the same approximation which neglects all the deriva- 
tives of the function fi. However, there the non-trivial 
"cross-horizon" correlations are contained entirely to the 
relationship between the lab coordinates x and t and the 
null co-ordinates x . Indeed also in [TO] these correla- 
tions across the horizon are demonstrated for an explic- 
itly 1+1 D black hole geometrjj^Jwhich, being conformally 
flat, possesses the standard logarithmic functional form 
for the correlator. Again the Hawking signal is contained 
exclusively in the relationship between the null and lab 
coordinate functions. 

The reason for this discrepancy is rooted in the rather 



7 Note the difference between the analyses of |10| and [8] [S]: 
whereas I10| works exclusively with the 1+1D conformally flat 
geometry the minimally coupled field the 1+1D action of IH1 19) 
arises from a dimensional reduction of a 3+ ID model. This di- 
mensional reduction is encoded in the 1+1D dynamics by the 
factor Q. which renders the field non-minimally coupled to the 
geometry. 



different causal structure of the black hole and cosmo- 
logical analogue spacetimes and in the inherently differ- 
ent natures of the associated particle creation processes. 
In the black hole spacetime the presence of the acoustic 
horizon is associated to an ergoregion which allows for 
negative energy states (w.r.t. an asymptotic observer) 
and these states in turn allow for the stationary Hawk- 
ing flux on the black hole static spacetime. In the case 
of the cosmological analogue geometry there is not such 
an ergoregion and particle production is permitted only 
as a result of the time dependence of the geometry. This 
time dependence is fully encoded in the conformal fac- 
tor describing the difference from Minkowski space of 
the analogue geometry. It is then clear that when one 
introduces an approximation which induces a conformal 
symmetry in the field equation, no particle production 
can take place. 



B. Spatially Homogeneous Q(Z) approximation 

As discussed above, the approximation in which we 
discard both temporal and spatial derivatives of f2 is too 
strong to capture the physics of particle creation in the 
1+1 dimensional massless case. However, we shall now 
argue that it is sufficient to keep the T derivatives while 
neglecting the Z derivatives of Vl to observe to particle 
creation signal. 

Firstly, since the relationship between Z and the co- 
moving lab coordinate z is given by integrating ( V. 1 ) 



(in the approximation which neglects radial and angular 
derivatives, i.e. / = 0) to 



Z(z) 



V2 



2/1: 



(V.12) 



we see that dz^l is naturally suppressed. One has 



Oz 
aZ 



(V.13) 



where the first factor dz/dZ is naturally small near the 
boundary (but still within the region in which the TF 
approximation remains valid) and the second factor d z £l 
is naturally small in the central region where the con- 
densate density is almost constant. Therefore, the factor 
dz^l is naturally small everywhere on the condensate. 

No such natural cancellation is available for <9yO since 
dt/dT = A 3 / 2 A ± , which is not naturally a small factor. 
Hence in what follows we will keep only T derivatives of 
n and discard the Z derivatives. 

Note that this approximation is not the same as as- 
suming to be constant in lab coordinates, the differ- 
ence being in the extra factor of dz/dZ which suppresses 
dz at the boundary of the condensate. An important 
feature of our analysis here is the inclusion of the effects 
of non-constant density and finite condensate size. 
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Phase-Phase Correlations 



Consider, then, (V.10) and assume ft is a function only 



of the variable T. Then (V.10) is written 



ft- 
7] 



dn 



(V.14) 



Let 77 be defined by the differential expression dT — f2 drj. 
Then (V.14) separates to 



d 2 ^ + n 2 (r,)k 2 k = 



(V.15) 



-ikZ 



is the fcth Fourier corn- 



where 4>k = (2tt)^ 1 J l dZ<pe 
ponent of the function <f> and L is the Z-svze of the con- 
densate. 

Since we work in the compact region L, inverse Fourier 
integrals will be replaced by discrete Fourier series in 
this section. Hence k in (V.15) is constrained to sat- 
isfy k = 2nir/L after choosing the boundary conditions 
4>(-L/2) = 0(L/2) = 0. 

Such an equation of motion can be solved exactly, as 
we did in the previous section whenever fi 2 is contained 
in the three parameter family of functions 



S 



tanh 



au a/, 770 G 



12 2 

(V.16) 

Our strategy here will be to make use of the available 
exact solutions and try to find a suitable element in S 
labelled by (a,, a f, 770) which most accurately reproduces 
the exact function Q associated with free expansion of the 
BEC. We will see that such a class of functions is indeed 
appropriate for an approximation to the free expansion 
case with the key benefit of a description in terms of 
asymptotic particle states. 

An alternative way to think about this approximation 
methodology is that one is studying the correlations in a 
BEC which really has this particular form for the func- 
tion n. As will be discussed later, there is a rather direct 
link between f2 and the scaling functions A. In fact from 



(V.6) we see that fl oc A so that fl € <S are mod- 



els for "finite expansion" condensates (or "finite contrac- 
tion" depending on the relative magnitudes of and 
a/). Such finite expansion (contraction) condensates are 
achievable in the lab with a suitable ramp-down or ramp- 



up trapping potential (see sec V C ) . Note that the finite 
expansion case qualitatively differs from the free expan- 
sion case only in the future asymptotically static region. 

Recall from Fig. [4] that the typical frequency of pro- 
duced particles during expansion is X(t)/X(t) which con- 
verges to zero at late times for a linear expansion A(t) cx t. 
For this reason we expect all the particle production in 
the infinite expansion case to occur only during the "ac- 
celerating" phase and that all of the interesting physics 
to have ceased by the time the scaling function becomes a 
linear function of time which occurs at late time. In this 
way we expect the finite expansion approximation to free 



expansion to capture the relevant physics neglecting only 
the irrelevant linear phase of expansion at late times. 

The strategy to understand particle production effects 
in the correlator is as follows. We choose the quantum 
state to be |in) (and since we work in the Heisenberg 
picture the system remains in this state for all time) and 
express our results in the rj — > 00 limit in terms of the 
number eigenstates of the asymptotic Hamiltonian which 
are the physical particle states in that region. 

Let us now move back to the 9 notation for the phase 
perturbation field. The field operator is written in 1+1D 
as the sum 



= XV*+&U 



(V.17) 



where the functions satisfy the equation of motion 
(V.15) and are proportional to plane waves in the past 



and linear combinations of plane waves in the future as 



fk(v — > 00, Z) cx sin kZ (a^e 



(V.18) 



The normalization for the functions fk is not determined 
by (V.15) but instead by the consistency between the 
commutation relations for with the commutation re- 
lations between the field operator 8 and p, the density 
perturbation field. To fix this normalization we will 
use the commutator [0(Z),dtO(Z')]. In 3+lD in the 
TF approximation with zero background velocity flow 
in the BEC one has p = —dtO/g as well as the funda- 
mental commutator [6(x), p(x')] = —ihS 3 (x 7 x') imply- 
ing [6(x),dt9(x')} — ihg5 3 (x, x'). Recall that we assume 
the perturbation field to be constant over the cross sec- 
tional area of the condensate (which was related to the 
existence of a UV cutoff at the healing length). Hence 
integrating this commutator we get for the 1+1D phase 
field an extra factor of the area from p( 3+1 ^Aj_ = p( 1+1 ^ 
and hence 

[0{x),d t 9{x')] = l ^-5{x,x'). (V.19) 

Further, expressing this commutator in Z coordinates we 
recall that the conjugate momentum is a density of weight 
one whereas the field 9 is a scalar so that one picks up a 
factor of the Jacobian under a coordinate transformation 



[9(Z),d t 9(Z')} = [9(x),d t 9(x')} x 



dZ 
dx 



Therefore 



[6{Z),d t 9{Z')] = ^f8{Z,Z') 



(V.20) 



(V.21) 



where c Si o is the initial sound velocity (which appears as 
a result of the use of c s instead of c s in the definition of 
the coordinate Z) and A/ is the final value of the scale 
factor A (which comes from the relationship between co- 
moving x and lab x). The dispersion relation is easily 
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computed to be uj 2 = £l 2 k 2 — > ajk 2 using the metric in 
(77, Z) coordinates. Hence we arrive at the expression 

ighXf 



S(Z,Z') = [0(Z),d t 6(Z')[ 



= £[fc.&l](/fcfc/*-/kftA) 

(V.22) 



fc,fc' 



which relates the commutators and the Wronskian for 
the mode functions allowing us to fix the normalization 
factor. 

For a Fock space (particle) representation we require 
[bk, b\., ] = 5kk> ■ Hence we arrive at the Wronskian con- 
straint on the mode functions 



E ^ kdt f'k' - fkdtfk') = 



ighXf 
c s nAj_ 



S(Z,Z'). (V.23) 



In order to take the laboratory time derivatives above we 
require the relationship between t and r\. We have 

dT 1 dt n 2 , 
*> = ir = n^ = Sg* (v - 24) 

where we have firstly transformed to the variable T using 
dT = fldi] and consequently to t using ( V. 1 ) , expressing 
the result in terms of Q alone. f2 and A are related by 
(V.6) as ft — £7o/\/A where flo is the spatial part of Q, 
numerically £Iq ~ 1-5 x 10 37 . Then 



dt 



f^Q drj 



tanh 



(V.25) 



Now, we wish to approximate the free expansion scale 
factor shown in Fig. [T] for approximately the first t = 
0.15s. This choice of the approximation region is moti- 
vated by two constraints: Firstly, after about t = 0.15s, a 
large proportion of the condensate no longer is described 
accurately by a TF approximation as can be seen in Fig. [3] 
; Secondly the natural shape of the exact free expansion 
scaling function is less accurately fit by the tanh func- 
tions we use here over longer time periods since the scal- 
ing function enters a linear regime at late times. 

Over this t = 0.15s time interval the scaling function 
A increases from Ai = 1 to approximately Xf — 25. Since 
fl = Oq/vA we see that the variable Q must decrease 
from Q to Clp/5 . Therefore we must have a/ = Qij/5. 
Hence by (V.25 1 we have the asymptotic relations be- 



tween the time variables 



t(v) = 



Xitt r], 
A/^o??, 



fort- 
fort 



—00 
00 



(V.26) 



where the function t(rj), obtained by integrating exactly 
the differential expression (V.25), varies smoothly be- 



tween these two constant asymptotes in the intermediate 



region. Taking the average gradient in the intermediate 
region for t(rf) we see that a t interval of 0.15s corresponds 
to an T) interval of approximately 12 x J7 ~ 8.2 x 10~ 40 
(in units of [r2 _1 ]s or [Energy] -[Time] 2 ). Therefore, we 
require the function O(ry) to vary over this time scale be- 
tween the initial and final values cii and af. That is we 
should choose 



.2 x 10~ 40 [Energy] • [Time] ; 



(V.27) 



With these choices one finds that ujkT] = kt/Xi 
For the un-normalized fk we have 



3/2 



E - fkdtfk) = 2 

k 

Scaling the functions fk as 

fk = 



E 



ik 

J/2 



sinfcZ sin/sZ'. 



xThg 



f 



Vol A±cik 



fk, 



(V.28) 



(V.29) 



where Vol is the axial Z-length of the condensate (which 
has the dimension of a time) , the Wronskian provides the 
correct numerical and functional form for the commuta- 
tor. 

Hence the correctly normalized mode functions are 
given by 



fk 



A/ /2 % sin(fcZ) - ikt /xr 



Vc*A ± Vol 



(V.30) 



Here we have transformed back to lab t coordinates in 
which we take the time derivatives. 

In 1+1D at late times the equal time Wightman func- 
tion is written in terms of otk and /3k in a very similar 
way to the 3+ ID example above in terms of the discrete 
(and inhomogencous) Fourier series 

00 

G e (i], Z, Z') =JVV- sin (k n Z) sin (k n Z') 
n 

n=l 

x (l + 2|/3 fe J 2 + 2$Hea fe7i /3£ e - 2lfe " Q ") 

(V.31) 

where the normalization factor N is given by 



TV = 



X f rig 

A±C s n%' 



(V.32) 



Taking symmetrically spaced points from the center of 
the condensate, the transient time dependent contribu- 
tion is given by 



1 = N 



OO 

E 

n=l 



1 . 


717T ^1 - 


2Z\ 






— sin 




sin 


nir 


n 











x 2D\eak n (3l e 



-2ik, 



2ZY 
(V.33) 
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FIG. 11. The bump of correlations propagating outwards 
from the center of the BEC cloud. The curves are for t = Os 
(red, solid), t = 0.3s (blue, dotted), t = 0.8s (green dashed) 
and t — 1.2s (black, dash dots) Here we plot the truncated 
sum up to the first 90 terms. The horizontal axis of this plot 
is the fraction of the half length L z /2 of the condensate after 
expansion. The vertical axis is the integrated contribution 
Ent in units of the normalization factor N. 



is easily computed as a truncated sum. 

With the choices made above for the parameters a*, 
a/ and 770 we find the entanglement structure shown in 
Fig. [TT] The time scales of the propagating bump are 
understood in terms of the crossing time for the mass- 
less modes: After the expansion the sound speed has 
decreased to one fifth its initial value c s /5 due to the di- 
lution of the BEC gas (A(t — > 00) = 25). The size of the 
condensate has increased by a factor of 25 so the crossing 
time increases by a factor of 125. Initially with a axial 
length of L z = 5 x 10~ 5 m and c s = 2 x 10 _3 ms _1 the 
crossing time is Crossing = 5 x 10~ 2 s which increases to the 
order of seconds after expansion. (Note again as in the 
3+ ID case, that the entanglement propagation speed is 
twice the sound speed since it represents pairs of phonons 
traveling in opposite directions). This is a slight overes- 
timation since the sound speed actually decreases from 
a maximum at the center of the condensate to zero at 
the edge. However we have neglected this slowing in the 
analysis. This order of magnitude estimate is confirmed 
in the plot of Fig. [TT] 

In Fig. [12] we show also the integrated contributions 
from the "1" and "|/3fc| 2 " terms in the integrand, in the 
same units as Fig. [TTJ We note that the magnitude of the 
bump structure is in general much larger than the back- 
ground ("1") and enhancement ("|/3fc| 2 ") contributions. 
This is expected since the produced particles propagate 
from small distances where there is a high intrinsic cor- 




0.012 



FIG. 12. The contributions to the phase Wightman function 
from the "1" and "|/3fc| 2 " terms in (V.31I truncated to the 



first 200 term of the series. Such a truncation serves as an ul- 
traviolet regulator on the correlation singularity at the origin 
of the "1" background contribution. No such UV divergence 
is present in the "|/3fc| 2 " contribution. 



relation to long distances where the natural background 
fluctuations are no longer correlated. 



2. Density- Density Correlations 

Experimentally, although there are proposals to di- 
rectly measure the phase correlations, the density corre- 
lations are those of major observational interest. Again, 
in the TF approximation the density field is written in 
terms of the phase as 



P 



(3+1 ^ = --d t 9 - -v z d z 



9 



9 



(V.34) 



The 1+1D density is related to the 3+1D density by 



( 3+1 ) = p( 1+1 'A± where A± is the cross sectional area 



16 



of the condensate (which we assume constant). Then the 
density correlator is derived according to 



G p (z, /) = \im t ^t>VG e \t, t', z, z) 
where we define the differential operator T> as 



V 



^ (d t + v z 8 z ) (d t > 
9 



(V.35) 



(V.36) 



In the 1+1D finite expansion case, T> reduces to sim- 
ply the time derivative since v z = after the expansion 
has finished. In Fig. [13] we plot the normalized density 
correlator G p (t,z,z') := G p (t, z, z')/(/0 (1+1) ) 2 for sym- 
metrically placed points z, z' about the middle of the 
condensate in units of N x A^/g 2 . We note again the 
propagating structure but of a slightly more complicated 
shape to the phase correlation bump. 



C. Finite expansion of an elongated BEC 

The availability of an asymptotically static regime in 
the future is a necessary ingredient for the application of 
the Bogoliubov formalism. As discussed above there are 
two ways to use the tools contained in that formalism 
to understand a system which does not possess such a 
future static regime such as the case of indefinite expan- 
sion of a released BEC. Either one considers the finite 
expansion as a theoretical approximation to the "true" 
unconstrained expansion case, as shown in Fig. [14] or 
alternatively one can imagine actually performing a fi- 
nite expansion experiment in the lab with a condensate 
which is at first released and consequently "caught" after 
a finite expansion. Above we have essentially been dis- 
cussing the first interpretation. Here we briefly discuss 
the second. 

In order for X(t) to follow the correct step up profile 
with a future non-dynamical region, one might expect 
a simple two level potential would suffice. However, the 



equations of motion (11.29) for the scaling function admit 



oscillating solutions so that in order to arrive at a static 
condensate in the asymptotic future a rather finely tuned 
potential is necessary. 

Again assuming the constancy of the radial scaling 
functions during the axial expansion we have for the trap- 
ping frequency 



wn(*) 



"jj(0) 
A»(f) 



Ajf) 
A(f)' 



(V.37) 



The necessary trapping frequency is obtained by insert- 
ing the desired expansion profile into this expression. 
One finds that in fact, the squared trapping frequency 
becomes negative for a short time near the end of the 
transition from expansion to staticity indicating a tem- 
porarily repulsive force is necessary to stop the BEC from 
contracting further. 



G(0.3,x)/N 




0.002 0.004 0.006 0.008 0.010 0.012 



G (0.8, x)/N 




0.002 0.004 0.006 0.008 0.010 0.012 



G(1.2,x)/N 



0.002 0.004 0.006 0.008 0.010 0.012 



FIG. 13. The density correlation function after expansion 
has taken place, for the finite expansion case, in units of N x 
A^/g 2 .. The curves from top to bottom correspond to t = 
0.3s, t = 0.8s and t = 1.2s after the onset of expansion (which 
we take, as described above, to last for 0.15s). 



In this way it would be possible to actually create in the 
lab a BEC system with particular time dependent trap- 
ping frequencies such that the scaling function A behaves 
in the way appropriate (as discussed above) to induce 
the kind of finite expansion or contraction analogue ge- 
ometry for perturbations. It is worth pointing out again 
at this point that we expect the general structure of a 
propagating bump of correlations resulting from a finite 
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that the condensate dynamics reduces to an effectively 
1 + 1-dimensional system and that the Wightman func- 



0.05 



0.10 



0.15 



FIG. 14. The exact scaling function A (red, solid) and its 
approximation by the tanh function Qq/Q (blue, dashed) over 
the interval t = 0.f5s as described in the text. 



scaling to not depend in any crucial way on the fine de- 
tails of the scaling. We expect the result to be insensitive 
for example to the fact that a precise reproduction of the 
tanh shape profile is practically very difficult and most 
probably would only be approximated in any real exper- 
iment . 



VI. CONCLUSION 

In this paper we have considered the correlation struc- 
ture of expanding BECs both in the idealized case of 
isotropic and homogeneous expansion as well as in the 
experimentally relevant case of a finite size, cigar-shaped, 
BEC anisotropically expanding. 

In the first case we have developed a general formalism 
for probing the effects of particle production on correla- 
tions as embodied in the equal-time Wightman function. 
Specifically, for any spatially translation invariant system 
we have shown that it makes sense to define the correla- 



tor as in ( III. 1 1 which in the case of particle production 



(using Bogoliubov coefficients) takes the form (IV.lll. 
This result is completely generic to particle production 
in any spatially translation invariant system, applying to 
standard cosmologies as much as the "analogue space- 
times" of central interest in the current article. Specifi- 
cally, the general cosmological result was shown to be di- 
rectly implementable in a 3+1D BEC system with a time 
dependent scattering length as well as a 3+ ID isotrop- 
ically expanding BEC, opening two alternative window 
on possible observation of this correlation signal. 

In the second case (anisotropic BEC) we have shown 



tion can be generically cast in the form ( V.31 1 for these 
systems. Using this result we have studied the correla- 
tion structure and shown how it carries the signature of 
analogue cosmological particle production possibly ob- 
servable in future experiments. In particular however, 
we have studied the density correlations in such BEC 
systems, which are of current experimental interest. 

The magnitude of the normalized correlations is of per- 
tinent experimental concern. The overall scale of the 
propagating bump is given in terms of the normalization 
factor for phase perturbation along with the additional 
factors which relate the 1+1D density to the phase field 
^|§^ ~ T^p^ijTa ~ 1CP 7 . Although this is a very small 
number, it should be noted that this magnitude is in fact 
much larger than the magnitude of the background corre- 
lation structure (present in flat spacetime also) . The rel- 
ative magnitude of the background to entanglement cor- 
relations is independent of the normalization of the mode 
function solutions but instead depends on the structure 
of the Bogoliubov coefficients which derive from the wave 
equation alone. 

The main purpose of this paper was to investigate the 
characteristic signature related to the cosmological par- 
ticle emission in an expanding (or more generally time- 
dependent) BEC, and analyze it without specific atten- 
tion to its actual measurability. However, one can of 
course imagine altering the parameters of the background 
BEC in such a way to improve the signal. For example, 
one parameter which can be altered by an order of mag- 
nitude is the number of atoms N which affects the value 
of po entering into the denominator of the normalized 
correlations. Another possibility of improving the sig- 
nal would be to alter the scattering length. We have 
discussed above the use of the analogue FRW geometry 
when one induces a time dependence in the scattering 
length. Here we have something different in mind. Since 
the background density depends on the scattering length 
as a 3 / 5 , by decreasing and holding the scattering length 
constant before expansion, one might be able to improve 
the signal. Another possibility is to choose a longer time 
scale on which to approximate the eternal expansion with 
a finite expansion; By stopping the expansion later, the 
correlations will be propagating on a background density 
which is lower which would improve the ratio in the nor- 
malized correlations. All these ways of altering po will 
also affect the values of 770 , on and a/ possibly giving a 
more pronounced signal. 

A separate analysis should regard the role of a non- 
vanishing initial temperature of the BEC (or more gen- 
erally the role of a non- vacuum initial state). In fact, as 
already pointed out in |23] , the net effect of such 
an initial temperature (i.e. a non "empty" initial state, 
but thermally populated) is to significantly increase the 
quantity of created excitations, therefore enhancing the 
signal without affecting the mechanism responsible for 
the particle production. We intend to make further in- 
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vestigations in these directions in future work. 
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